The continuum limit of quark number susceptibilities 
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I. INTRODUCTION 
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We report the continuum limit of quark number susceptibilities in quenched QCD. Deviations 
from ideal gas behaviour at temperature T increase as the lattice spacing is decreased from T/4 
to T/6, but a further decrease seems to have very little effect. The measured susceptibilities are 
20% lower than the ideal gas values, and also 10% below the hard thermal loop (HTL) results. The 
off-diagonal susceptibility is several orders of magnitude smaller than the HTL results. We verify 
a strong correlation between the lowest screening mass and the susceptibility. We also show that 
the quark number susceptibilities give a reasonable account of the Wroblewski parameter, which 
measures the strangeness yield in a heavy-ion collision. 
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With the RHIC now in its second year of running, it is important to pin down lattice predictions for the high 
temperature phase of QCD. In this context, fully non-perturbative measurements of quark number susceptibilities 
jl], || ||, ^| are important for four reasons. Firstly, they are directly related to experimental measurements of event-to- 
\Q | event fluctuations in particle production ||. Secondly, earlier results ^| showed a strong jump in the susceptibility 
■ across the phase transition, but indicated a significant departure from weak-coupling behaviour, and it is important 
to check whether this persists into the continuum. Third, resummed perturbative computations of this quantity have 
now become available M , thus making it possible to accurately test the importance of non-perturbative contributions 
to this quantity. Finally, with continuum extrapolated results in hand, one can address the question of whether the 
strangeness production seen in heavy-ion collisions can be quantitatively explained as a signal of the quark-gluon 
plasma. 

We have recently presented systematic results for quark number susceptibilities in quenched QCD || as well as for 
QCD with two flavours of light dynamical quarks [Q. In these studies a large range of temperatures, T, was covered 
O i| at a series of different quark masses, to. These computations were done at a fixed cutoff with the lattice spacing, 
a = T /4 while keeping the finite volume effects under control so that the thermodynamic limit could be taken reliably. 
' There was a 3-5% difference between the quenched and dynamical computations. In this paper we examine the cutoff 
^ . dependence of the quenched results and report their continuum (zero lattice spacing) limit. 
The partition function of QCD with three quark flavours is 
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where the temperature T determines the size of the Euclidean time direction, S(T) is the gluonic part of the action and 
the determinants of the Dirac matrices, M, contain as parameters the quark masses, to/, and the chemical potentials 
/// for each of the flavours / = u,d, s. We also define the chemical potentials /io = fi u + fJ-d + Ms; M3 = M« ~ Md an d 
Ms = Mm + Md ~ 2/i s , which correspond to the diagonal flavour SU(3) generators. Note that fio is the usual baryon 
chemical potential and ^3 is an isovector chemical potential. 
Quark number densities are defined as 

T dlogZ 

ni{T,^ u ,fi d ,^ s ) = —— , (2) 

V d[ii 



'Electronic address: ;avai@tifr.rcs.in 
^Electronic address: sgupta@tifr.res.in 



2 



and the susceptibilities as 

T OHogZ 

Xy (T, M , Ms ) = -^-, (3) 

where V denotes the spatial volume. The subscripts i and j are either of the index sets / or a — 0, 3 and 8. 
We lighten the notation by writing the diagonal susceptibilities xu as Xi- We determine the susceptibilities at zero 
chemical potential, fif = 0. In this limit, all rii(T) — 0. Since we work with m u = rrid = m < m s , we also have 
X03 = 0. 

Flavour off-diagonal susceptibilities such as 

Xud = (^r) (trM-'Mi trM^M'a) (4) 

are given entirely in terms of the expectation values of disconnected loops. Since m u = md, we obtain \us — Xds with 
each defined by an obvious generalization of the formula above. Of the flavour diagonal susceptibilities we shall use 



(f ) [(( trM *" lM «) 2 ) + ( tr (M^K' - M^MiM^Ml))] . (5) 



Xu — Xd are given by a generalization of this formula. Numerically, the simplest quantity to evaluate is the diagonal 
iso- vector susceptibility 

X3 = \ (tr (M^M: - M^M' U M^M' U ) ) . (6) 

Two more susceptibilities are of interest. These are the baryon number and electric charge susceptibilities, 

Xo = ^(4x3 + Xs+4x„ d + 4xn S ) and \ q = ^ (10X3 + Xs + Xud - ^Xus) ■ (7) 

Note that \o is the baryon number susceptibility for three flavours of quarks. As a result, this expression differs from 
the iso-singlet quark number susceptibility for two flavours, defined in [lj], both in overall normalization and by terms 
containing strangeness. In our numerical work we have chosen to use staggered quarks. Hence, to normalise to one 
flavour of continuum quarks and compensate for fermion doubling on the lattice, we have to multiply each of the 
traces in eqs. (||-|6|) by a factor of 1/4. 

The quark mass appears in two logically distinct places — first in the operators which define the susceptibilities 
in eqs. (|^"|^), and secondly in the determinant in the partition function of eq. (Q) which defines the weight for the 
averaging of these operators. The first defines a valence quark mass and the second the sea quark mass. In this paper 
we adopt the quenched approximation, whereby the determinants in eq. (Q) are set equal to unity. However, we shall 
vary the valence quark masses. As shown earlier, unquenching the sea quarks changes results by 3-5% || 0]. 

Note that all the flavour off-diagonal susceptibilities are exactly zero in an ideal gas. We shall denote by x~fft the 
ideal gas value for \3- On an N t x lattice with N c colours and lattice spacing a, we find 



° 2 Xfft = ^J^Y,{ D 2 s™ 2 po cos 2 Po + D 1 (sin 2 p - ~ J j , 



(8) 



where the spectrum of momenta is po = (2n/Nt)(n + 1/2) with < n < N t and p t = 2ttti/N s with < n < N s , and 
D = (ma) 2 + sin 2 p^. For a given m/T c and T/T c , the value of the quark mass in lattice units is 



N t ) \TJ (t 



(9) 



The formula in eq. (g) is normalised such that there is exactly one copy of each flavour of quarks. The values of \ s , 
Xo and Xq f° r an ideal gas can then be simply obtained as linear combinations of Xfft f° r different valence quark 
masses. 
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FIG. 1: The temperature for any value of Nt for a given coupling /3. Crosses denote data for N t — 4 lattices and pluses for 
Nt = 8 lattices. The scale has been set by a corrected QCD 2-loop formula The difference between the two sets of data is 
a measure of scale breaking by residual power-law corrections and is negligible within statistical errors. 



II. LATTICE RESULTS 



We have previously reported measurements of \3 m quenched QCD with 4 time slices (Nt = 4) at temperatures 
of 1.1T C , 1.25T C , 1.5T C , 2T C and 3T C Q. We have extended these measurements to higher temperatures through two 
simulations on 4 x 20 3 lattices at (5 = 6.6 and 6.7. These correspond to temperatures of 4.13T C and 4.7T C . All these 
computations are performed at fixed bare quark masses of m/T c = 0.03, 0.3, 0.5, 0.75 and 1 (the data for m/T c = 0.5 
are new). 

The continuum limit has been taken by going to several smaller lattice spacings. One set of computations is 
performed with a lattice spacing a = 1/6T which is 33% smaller than the lattice spacing used earlier. We have 
taken data on 6 x 20 3 lattices at /3 ~ 6.0625, 6.3384 and 6.7 corresponding to T/T c = 1.333, 2 and 3.133, and on a 
6 x 24 3 lattice at (3 = 6.8 corresponding to T/T c = 3.532. One further computation was made on an 8 x 18 3 lattice at 
(3 = 6.55. This corresponds to T/T c = 2 with lattice spacing a = 1/8T. The quark mass m/T c is kept independent of 
a and T, and hence ma decreases with increasing T at fixed N t or with increasing Nt at fixed T according to eq. (|^). 

The lattice scale has been set using the plaquette measurements of || and the analysis performed in ||. The 
conversion of (3 to T using data for Nt = 4 and 8 give the two curves in Figure [|. The difference between the two 
curves is within the statistical uncertainty of 5% on the scale assignment. There is the same degree of statistical 
uncertainty in the bare quark masses and in all other scales. We shall not display this inherent scale uncertainty in 
the measurements, but it should be kept in mind. 

The simulations have been performed with a Cabbibo-Marinari pseudo-heatbath technique with 3 SU(2) subgroups 
updated on each hit. An initial 1000 sweeps have been discarded for thermalisation. On the N t < 6 lattices we 
have used 80 configurations separated by 1000 sweeps for the measurements. On the N t = 8 lattice we have used 55 
configurations separated by 500 sweeps. 

We have previously checked that the spatial lattice size, L — N s a, is quite irrelevant to the values of x measured as 
long as it is sufficiently large compared to the inverse pion mass |2|, 0|. In fact, at all the couplings and quark masses 
we have used, m^L > 5. The most important other constraint on L comes from the requirement that it should be 
large enough to prevent spatial deconfinement 0. This is ensured by taking N s /N t > T/T c in all our simulations. 
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FIG. 2: Xs/Xfft at T = 2T C shown as a function of a 
evaluated on lattices with N t = 6. 



oc 1/Nt for m/Tc = 0.03. The ratio reaches its continuum limit when 



The traces are evaluated by the usual stochastic technique, 



Tr 



1 N 

2TV ^ 1 



(10) 



where Ri are a set of TV uncorrelated vectors with components drawn independently from a Gaussian ensemble with 
unit variance. Each vector has three colour components at each site of the lattice. Since we use a half lattice version 
of the Dirac operator for staggered Fermions, the number of components of each vector Ri is one and a half times 
the number of lattice points. (Tr A) 2 is obtained by dividing the set of Ri into disjoint blocks, constructing Tr A in 
each block, taking all possible distinct pairs of such estimates, multiplying them and then averaging over the pairs. 
For TV t — 4 it is possible to get accurate results with TV « 10, although we have chosen to use TV" = 80 for each gauge 
configuration 

It is easy to check that the variance of the estimator above is proportional to Tr(A 2 ). Since the diagonal element of 
the Fermion matrix is proportional to the quark mass ma, whereas the off-diagonal elements are bounded by unity, 
with increasing T or decreasing a, Tr A decreases linearly, while its variance remains constant. Hence, for sufficiently 
small ma the number of vectors has to increase quadratically with ma. We found that TV„ = 100 was sufficient for 
N t — 8 at T = 2T C , although at larger N t or T significantly larger values of TV,, are needed. A further numerical 
problem arises from the cancellation of the two matrix elements which give %3 (see eq. |6|) . While increasing TV t and 
keeping the physical size of the lattice fixed, each of these terms increases quadratically but mutually cancel to give a 
number which decreases quadratically with TV t . At fixed word length, this implies a reduction in accuracy by a factor 

<>r V. 

In Figure |2| we show how the continuum limit of \3 is reached. There is an 8% decrease in \3 1 Xfft m S om g from 
TV t = 4 to TV t = 6 lattices. However, within measurement errors, the continuum limit is reached already at TV t = 6. 
This implies a cancellation of the TV t -dependence of %3 and Xfft- ^ ne reason for the atypical behaviour for TV t = 4 
is not in the numerator, which is non-perturbative, but in the denominator, which is the ideal gas. 

For TVt = 4, the spectrum of po consists of {±7r/4, ±37r/4}. As a result, sin 2 po = 1/2 and the second term in the 
sum on the right of eq. (JsJ) vanishes identically. This accident occurs only for TV t = 4, and for all other values of TV t 
this term contributes a non-vanishing value. It is easy to see that in the limit TV t — > oo this term contributes to the 
leading part of x%ft- This is the reason the values for X3/Xfft f° r TVt = 4 lie well away from the continuum limit. 
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FIG. 3: The ratio Xz/Xfft i n quenched QCD, shown as a function of T for the values of m/T c indicated. Pluses denote data 
for lattice spacing a — 1/4T and circles for a = 1/6T. The lines are cubic spline fits to the data. 



T/T c 


X'a/XfFT 




d /T 2 


Xs/Xfft 


Xo/xfft 


X.q/X.FFT 


1.333 


0.70 (3) 


(2 ±4) 


X 


10" 6 


0.455 (5) 


0.66 (2) 


0.68 (3) 


2.000 


0.761 (3) 


(-1±2) 


X 


lO" 7 


0.656 (3) 


0.740 (3) 


0.751 (3) 


3.133 


0.801 (2) 


(-1±1) 


X 


io- 7 


0.766 (2) 


0.794 (2) 


0.798 (2) 


3.532 


0.807 (2) 


(1±2) 


X 


lO" 7 


0.789 (2) 


0.803 (2) 


0.805 (2) 



TABLE I: Results for the continuum limit of quark number susceptibility in quenched QCD. We have taken m a /T c 
appropriate to full QCD. Our N t = 4 results indicate that at higher T the ratios x/Xfft remain flat. 
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Our complete results for the dependence of Xz/Xfft on ^ jT c at several different quark masses for N t — 4 and 6 
are shown in Figure ||. In view of the data shown in Figure our results for N t = 6 are the continuum limit for the 
quenched theory: these are summarised in Table |[ Note that even at temperature as high as 3.5T C the result is far 
below the ideal gas limit. Closer to T c there is an even more dramatic fall in the value of X3/Xfft- We have also 
measured the off-diagonal susceptibility Xud- As already seen for N t — 4 this is consistent with zero. In Table | we 
show that the continuum limit of Xud also vanishes within errors of 10 -7 T 2 . It is therefore much smaller than the 
3- loop perturbative result given in S . 

Note that the continuum limit of our measurements lie 20-30% below the ideal gas result. They also differ signifi- 
cantly from HTL as well as the skeleton graph resummed results of M. At T = 3T C , HTL predicts Xz/Xfft = 0.90- 
0.94 on varying the scale of a s between irT and 47rT, and the resummed computation gives x^/Xfft = 0.95-0.97. 
Our measurement shows Xz/Xfft = 0-80 at 3T C (see Table |). A second computation of the HTL result is available for 
Nf = 2 JlCfl . Since this agrees with the HTL result of Q for Nf — 2, we believe the HTL result for Xz/Xfft 1S under 
control. There is thus a genuine discrepancy between these lattice results and existing perturbative computations. 

The quark number susceptibilities are closely related to the screening correlator of a one-link separated quark 
bilinear operator which corresponds to a p meson at zero temperature. The T > transfer matrix of the problem 
mixes this with a quark bilinear that corresponds to the n [yj. We have earlier shown evidence for Nt = 4 that X3 
is closely related to the only known non-perturbative quantity among the screening correlators, that is the screening 
mass Ms coming from the two degenerate correlators that descend from the it and a operators of the zero temperature 
theory M. In view of the strongly non-perturbative character of X3: as revealed by the comparison with HTL and 
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FIG. 4: Correlation between the screening mass in the scalar/pseudo-scalar sector, Ms, and the quark number susceptibility, 
X'i at 2T C . The measurements have been made on a 8 x 18 3 lattice at a series of quark masses. The line is an exponential fit 
to the data. 



resummation explained above, examining this correlation in the continuum becomes more significant. In Figure |4| we 
show data obtained on 8 x 18 3 lattices which indicate that this correlation survives into the continuum. 



III. STRANGENESS 

Measurements of hadron yields in relativistic heavy-ion collisions have been analysed extensively [[[2], and the 
observed enhancement of strangeness has been claimed to be evidence for quark-gluon plasma formation. The chemical 
composition is observed at the hadron freezeout temperature, which has been found to be very close to the value of 
T c in full QCD. It has been argued that this chemical composition cannot arise due to hadronic rescattering It 
would be interesting if it could be directly determined whether or not the composition is characteristic of the plasma 
for T > T c . We argue below that this question can be answered by a lattice QCD measurement such as we have 
presented above. 

The quantity which seems to be directly connected to lattice measurements is the Wroblewski parameter |Iq| , X s , 
which measures the ratio of newly created primary strange to light quarks — 

A,= 2 ^- x . (11) 
(uu + da) 

A fluctuation-dissipation theorem relates the rate of production of q uark pairs by the plasma in equilibrium to the 
imaginary part of the full complex quark number susceptibilities . If the inverse of the characteristic time scales 
of the QCD plasma are not close to the typical energy scales for the production of strange and lighter quarks, then 
the ratio of their production rates is just the ratio of the static susceptibilities that we have measured. Then, if the 
observed chemical composition is created in equilibrium, we should have 

A s = « (12) 

Xu + Xd Xu 

where the last equality holds in the limit of equal u and d quark masses. These susceptibilities have to be evaluated 
on the lattice at the temperature and chemical potential, /xq, relevant to the collision. However, it has been shown 
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FIG. 5: The Wroblewski parameter, A s , obtained from various experimental measurements, compared with the value obtained 
from the lattice measurements reported here. 



recently |L7| that X s is insensitive to /io for SpS and RHIC energies. Furthermore a lattice measurement shows that the 
screening mass, Ms, does not change rapidly with increasing chemical potential fisj ]. Then, in view of the correlation 
between Ms and x shown in Figure |J, we do not expect rapid changes with chemical potential in the ratio above. 
Both these arguments show that as a first approximation one can take the ratio in eq. ( |l2| ) at zero chemical potential. 

From our measurements of x reported here, we can form the lattice "prediction" of A s created in the plasma at T c . 
This has been done by taking our measurements at T/T c = 4/3 and 2, and extrapolating them linearly to T c . For Xu 
we use data taken for m/T c = 0.03, which is light enough for this purpose. Since scaled quantities feel the smallest 
correction for unquenching 0, we take Xs to correspond to m/T c = 1, which is appropriate to full QCD. In Figure |5| 
the resultant prediction is compared with the data collected in [[l3[ |l7| . It can be seen that there is a fair agreement 
with the data. However, this is subject to several assumptions — 

1. The foremost assumption is that the characteristic time scales of the plasma are not close to the inverse energy 
scale of the production processes. It is not possible to test this assumption in an Euclidean computation. 
However, it has been suggested that these characteristic time scales could be observed in real or virtual photon 
production as spikes in the spectrum |]l9| . If these are seen and they lie in the energy range for the production, 
the assumption would be falsified. 

2. Another important assumption is that chemical equilibration takes place in the plasma. If chemical equilibrium 
is not achieved but the light degrees of freedom achieve energy equipartition, then the effective temperature 
could be higher. 

3. We have assumed slow variation of the ratio Xs/Xu with chemical potential. While this seems to be a reasonable 
estimate in view of the results of |[7| and our own observations about the relation between the lowest 
screening mass and the susceptibilities, one should keep in mind possible changes in the ratio due to violation 
of this assumption. 

4. Our continuum results are obtained in quenched QCD, whereas the data corresponds to full QCD. As mentioned 
already, we have taken m/T c — 1, as appropriate to full QCD, in order to correct for this. However, an additional 
5-10% shift in our results could arise when unquenching This is not shown in Figure |^. 
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5. Our measurements are made away from T c and extrapolated down to this temperature. In principle this can 
be corrected by a computation made directly at T c . However, due to the enormous cost of extracting physical 
quantities at T c , we postpone this to the time when we attempt a full QCD simulation. 

It might seem surprising that X s corresponds to the chemical composition at T c . However, a moment's thought shows 
that this is entirely in accord with the assumption of chemical equilibrium. Even if the system is thermalized at a 
much higher temperature, it remains in equilibrium until it reaches T c . Departure from equilibrium as the system 
cools through the transition is unlikely, and if it occurs, should have a visible signature such as the formation of a 
disoriented chiral condensate |2C[| . In the absence of such phenomena, if X s does not correspond to the chemistry 
at T c , then the second assumption above would be demonstrated to be false. Also, with increasing ion energy the 
system moves closer to zero baryon number, as seen at RHIC. This makes our third assumption better at the RHIC 
and LHC. The fourth assumption can be removed by a computation with full QCD. Such a computation is planned, 
where also the final assumption can be removed, and will be reported elsewhere. 



IV. SUMMARY 



In summary, we have presented new and precise results on quark number susceptibilities over a wide range of 
temperatures and quark masses in the high temperature phase of QCD. The main aim of this study was to obtain the 
continuum extrapolation of these susceptibilities. As shown in Figure |2| while there is a significant change in going 
from N t — 4 to N t — 6 lattices, there is no statistically significant change in the susceptibilities beyond this. The 
results obtained for N t — 6 lattices can then be taken as the continuum limit. This limit is shown in Figure |^ and 
Table |. 

There is a strong discrepancy between the continuum extrapolated lattice results and HTL computations for these 
quantities — varying between about 12% at the highest T to even larger at smaller T. The off-diagonal susceptibility, 
Xud-, is zero, in contrast to the HTL results. There is a somewhat larger discrepancy between the continuum extrapo- 
lated lattice results and the skeleton graph resummed results. The conjecture that there is a strong non-perturbative 
component to the quark number susceptibilities is supported by an observed strong correlation between the smallest 
quark bilinear screening mass, M$, and the susceptibility %3 (Figure [I]). 

It is interesting to note that the continuum extrapolated results for the strange and light quark susceptibilities can 
be used to give a surprisingly good description of the chemical composition of hadrons at freezeout in SpS and RHIC 
experiments (Figure 0). This has to be treated as a preliminary estimate due to the many caveats which we have 
listed, and some of which we plan to check in future lattice measurements. 
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